######################################################################
#Replication file for "The Impact of China's AIIB on the World Bank"
#Authors: Jing Qian, James Vreeland, Jianzhi Zhao
#Codes to replicate results in Table A.1, column (1) - column (4) 
######################################################################
rm(list = ls())


#Other packages used
library(tidyverse)
library(lmtest)
library(sandwich)
library(MASS)
library(stargazer)

#Setwd
setwd("C:/Users/qianj/Dropbox (Princeton)/AIIB_WB_Replication")

#Load custom functions
#Note: Given the number of tests performed in the appendix, 
      #several custom functions are created to implement estimation methods in a more efficient way.
      #These functions merely create a wrap-up of existing functions, like bpCausal and gsynth,
      #in order to avoid typing parameters repeatedly.

source("code/custom_functions.R")



###################################
#Table A.1, column (1) - column (4)
###################################

#-----------------------
#(1) Load data
#-----------------------
load("data/data_main.RData")

#-----------------------
#(2) Setup
#-----------------------
#Covariates
covs.all = c("gdppc_log_lag",
             "population_log_lag",
             "election_either_lag",
             "fdi_gdp_lag",
             "debt_gni_lag",
             "oda_gni_lag",
             "polity2_lag",
             "unsc_lag",
             "IdealPointDistance_lag")

#DV
dv.use = "hard_project_count_all"

#FE
fe.use = c("factor(iso2c)",
           "factor(year)")

#Post 2016 dummy
df.main$post_2016 = ifelse(df.main$year >= 2016, 1, 0)

#DiD Vars
var.did = c("aiib_founder_2016",
            "aiib_founder",
            "post_2016")

#Data to use
df.use = df.main

#-----------------------------
#(3) Column (1) Difference-in-Differences
#-----------------------------
#Create equation
equa = make.equation(dv.use,
                     c(var.did, covs.all))

#Estimate
fit.tab.a1.1 = lm(equa, data = df.use)
se.tab.a1.1 = fit.tab.a1.1 %>%
  coeftest(., 
           vcov = vcovCL(., 
                       cluster = df.use$iso2c))

#-----------------------------
#(4) Column (2) Negative Binomial
#-----------------------------
#Create equation
equa = make.equation(dv.use,
                     c("aiib_founder_2016", covs.all, fe.use))

#Estimate
fit.tab.a1.2 = glm.nb(equa, data = df.use)
se.tab.a1.2 = fit.tab.a1.2 %>%
  coeftest(., 
           vcov = vcovCL(., 
                         cluster = df.use$iso2c))

#-----------------------------
#(5) Column (3) Poisson
#-----------------------------
#Create equation
equa = make.equation(dv.use,
                     c("aiib_founder_2016", covs.all, fe.use))

#Estimate
fit.tab.a1.3 = glm(equa, data = df.use, family = "poisson")
se.tab.a1.3 = fit.tab.a1.3 %>%
  coeftest(., 
           vcov = vcovCL(., 
                         cluster = df.use$iso2c))

#------------------------------
#(6) Column (4), Two-way FEs
#------------------------------
#Create equation
equa = make.equation(dv.use,
                     c("aiib_founder_2016", covs.all, fe.use))

#Estimate
fit.tab.a1.4 = lm(equa, data = df.use)
se.tab.a1.4 = fit.tab.a1.4 %>%
  coeftest(.,
           vcov = vcovCL(., 
                         cluster = df.use$iso2c))

#------------------------------
#(7) Combined results for Table A.1, column (1) - column (4)
#------------------------------
stargazer(list(fit.tab.a1.1,
               fit.tab.a1.2,
               fit.tab.a1.3,
               fit.tab.a1.4),
          se = list(se.tab.a1.1[,2],
                    se.tab.a1.2[,2],
                    se.tab.a1.3[,2],
                    se.tab.a1.4[,2]),
          omit = c("iso2c", "year", "Constant"),
          type = "text")

